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Abstract 

In control and engineering community, models generally contain 
a number of parameters which are unknown or roughly known. A 
complete knowledge of these parameters is critical to describe and 
analyze the dynamics of the system. This paper develops a novel ap¬ 
proach of estimating unknown time-varying parameters of nonlinear 
systems based on real-time nonlinear receding horizon control strat¬ 
egy. For this purpose, the problem of estimation is converted into 
optimization problem which can be solved by backwards sweep Ric- 
cati method. Methodology is implemented on two nonlinear examples 
without providing any reference of the parameters. Results success¬ 
fully demonstrate the fact that proposed algorithm is able to estimate 
unknown time-varying parameters of nonlinear systems effectively. 


1 Introduction 

Parameter estimation is a widely used process in modeling of dynamic sys¬ 
tems, where systems exhibit unknown but desired dynamics. In real world 
applications, many dynamical systems, including systems which belong to 
the fields of biology, chemistry, physics, engineering and many more [1-8], 
may (and mostly do) have partial or all unknown parameters. Knowledge 


1 


about the parameters of the underlying dynamics is the prerequisite to an¬ 
alyze, control, and predict their characteristic behaviors. Hence, this topic 
has drawn great attention in various areas due to its theoretical and practical 
significance. 

While estimating unknown parameters in dynamic systems could be a 
daunting task, study of parameter estimation, especially in nonlinear systems 
that exhibit time-varying behavior, remains an open held of research. In this 
area, an important effort is the application of adaptive control and synchro¬ 
nization methodologies by designing an adaptive controller to estimate the 
uncertainties as well as minimize the synchronization error [9-14], However, 
adaptive control based methodologies highly depend on the assumption that 
parameters are constant or slowly time-varying. In practice, the time-varying 
feature of system parameters can have significant impact on the performance 
of the underlying systems [15-18], and canot be underestimated. 

On the other hand, especially with the improvements in microprocessors 
and increased computational capabilities, Receding Horizon Control (RHC) 
based methodologies gained significant momentum in the last two decades 
[19-24,26]. RHC is a formulation which is widely used to obtain an optimal 
feedback control law by minimizing the desired performance index for a given 
finite horizon. The performance index of a receding-horizon control problem 
has a moving initial time and a moving terminal time, and the time interval 
of the performance index is finite. Since the time interval of the performance 
index is finite, the optimal feedback law can be determined even for a system 
that is not stabilizable. One advantage of receding horizon optimal control 
formulation is that it can deal with a broader class of control objectives other 
than asymptotic stabilization. Moreover, nonlinear receding horizon control 
(NRHC) has made a significant impact on industrial controls and is being 
increasingly applied in process controls [23,24,26]. Various advantages are 
known for NRHC, including the ability to handle time-varying and nonlinear 
systems, input/output constraints and plant/parameter related uncertainty. 

There are many existing valuable works related to the use of different 
methodologies in nonlinear systems for parameter estimation procedures. 
One of the conducted studies concentrates on parameter identification us¬ 
ing differential evolution algorithms [2], while another study provides a bio¬ 
geography based parameter estimation for nonlinear systems [3]. Other im¬ 
portant studies include Kalman filtering based parameter estimation [5], 
adaptive estimation with invariant manifolds [7], estimation through iden¬ 
tical or non-identical structures [14], chaotic ant-swarm models [27,28], time 


2 


series approach [29,30], via particle swarm optimization [31], recursive iden¬ 
tification via incremental estimation [32], 

The main goal in this approach remains as the parameter estimation of 
nonlinear time varying systems, but the fundamental difference of this study 
with respect to the existing literature is that in this paper it is desired to 
demonstrate the connection (and applicability) of real-time nonlinear reced¬ 
ing horizon control algorithms as an effective parameter estimation proce¬ 
dure. Here, we provide a framework which is able to cope with nonlinear 
time-varying systems. For this purpose, the estimation procedure is reduced 
to a family of finite horizon optimization problems. In addition, a non¬ 
iterative optimization algorithm is employed to avoid high computational 
complexity and large data storage. With this formulation, we also provide 
the closed-loop asymptotic stability guarantees of the synchronization error 
system. Obtained results from two example problems demonstrate the ap¬ 
plicability of the real-time NRHC as a parameter estimation routine on the 
nonlinear systems with unknown time-varying parameters. This approach 
also serves as a minimization routine the synchronization errors, as well. 

The paper is organized as the following: Problem formulation is defined 
in Section-2. In Section-3, backward-sweep method is iterated and stability 
proof is provided. In Section-4, we discuss the findings of such approach, and 
with the conclusions section, we finalize the paper. 


2 Problem formulation 

For this study, we consider the following nonlinear system with unknown 
time-varying parameters: 

x = Ax + f{x) + D(x)Q(t), (1) 

where x G R n is the state vector, A G R nxn is the the linear coefficient matrix 
and f(x ) : R n —>• R n is the nonlinear part of system in (1). D(x) : R n —> R nx P 
is a known function vector and 0 G R p denotes the unknown time-varying 
parameters. 

In order to estimate the unknown parameters, the corresponding con¬ 
trolled system is given by 

y = By + f(y) + D(y)Q(t), (2) 
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where y G R n is the state vector and 0 represents the estimated parameter. 
Here, function /(■) and D(-) satisfy the global Lipschitz condition. Therefore 
there exist positive constants a and f3 such that 

11/(2/) - /0)ll <Pi\\y-x\\, 

\\D(y) - D(x)\\ <P 2 \\y-x\\. 

The system in Eq.(l) is considered as the drive system and system in 
Eq.(2) is considered as the response system. The synchronization error be¬ 
tween the drive and the response system becomes an important part of this 
analysis, and in this study it is defined as 

e(t) = y{t) -x(t). 

With this definition, it is possible to define the synchronization error system 
as 

e(t) = Ae(t) + B(f(y(t)) - f(x(t ))) + D(y)Q(t) - D(x)Q(t) (3) 
where the estimation error is denoted by 

0(t) = ©(*) - 0(t). 

From this point on, for the estimation of unknown time-varying parame¬ 
ters, we propose and establish a finite horizon cost function which is associ¬ 
ated with the synchronization error and the estimated parameter, as: 

1 f t+T 

J = - j [e T Qe + 0 T i?0]dr, (4) 

where weighting matrices Q > 0 and R > 0. With the construction of the 
synchronization cost function, as given in Eq.(4), the estimation problem is 
converted into a parameter optimization procedure, where the unknown pa¬ 
rameter^) can be estimated through this process. For this purpose, we utilize 
the powerful nature of real-time nonlinear receding horizon control algorithm 
to minimize associated synchronization the cost function. In this context, 
the performance index evaluates the performance from the present time t to 
the finite future t + T, and is minimized for each time segment t starting from 
y(t). With this structure, it is possible to convert the present receding hori¬ 
zon control problem into a family of finite horizon optimal control problems 
on the artificial r axis parameterized by time t. 
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It is well known from literature that first order necessary conditions of 
optimality are obtained from the two-point boundary value problem (TP- 
BVP) [19], as 


y* = H L y* =y(t), = 

X*(T,t) = 0, Hq = 0 
where H is the Hamiltonian and is defined as 
H = L + X * T y 

= \( eT Q e + ® T R9) + A * T [Ay + f(y) + D(y)Q(t)]. ^ ^ 

In Eqs.(5)-(6), ( )* denotes a variable in the optimal control problem so as 
to distinguish it from its correspondence in the original problem. 

Using this formulation, the estimated parameter can be calculated as 

©(f) = arg{H$[y(t), X(t),Q(t),x(t)] =0} (7) 

In this methodology, the TPBVP is to be regarded as a nonlinear equation 
with respect to the costate at r = 0 as 


F(\{t),y(t),T,t) = \*(T,t) = 0. (8) 


Since the nonlinear equation F(\(t),y(t),T,t) has to be satisfied at any 
time t, = 0 holds along the trajectory of the closed-loop system of the 
receding horizon control. If T is a smooth function of time t, the solution 
of F(X(t),y(t),T(t),t) can be tracked with respect to time. In this formula¬ 
tion, the ordinary differential equation of A (t) can be solved numerically, in 
real-time, without any need of an iterative optimization routine. However, 
numerical errors associated with the solution may accumulate as the inte¬ 
gration proceeds in practice, and therefore some correction techniques are 
required to correct such errors in the solution. To address this problem, a 
stabilized continuation method [33-36] is used. According to this method, it 
is possible to rewrite the statement as 


d F 
d t 


-A S F , 


(9) 
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where A s denotes a stable matrix to make the solution converge to zero 
exponentially, such that 

F = F 0 e~ Aat ->■ 0 (t ->• oo), (10) 

where Fq denotes the initial value of F. And thus the numerical errors is 
attenuated through the integration process. 


3 Backward-sweep Method: 


In order to integrate the differential equation of A(£) in real time, the partial 
differentation of Eqs. (5) with respect to time t and r, converts the problem 
to the following linear differential equation: 


d_ 

where G — f y — L = K = H yy - Since 

Xt and x T are canceled in Eq.(ll), data storage is reduced. 

The derivative of the nonlinear function F with respect to time is given 
by 

F = + K(T,t)F (1 2) 

In order to reduce the computational cost without resourcing any approx¬ 
imation technique, the backward-sweep method is implemented. This pro¬ 
vides the relationship between the costate and other variables is expressed 
as the followings: 


y* t - y* 


i 

i 


y*t - y* 

K-K 


i 

1 

>! 

1 

O 

4 


a * - a; 


( 11 ) 


K ~ K = s ( t, t ) (y* t -y*) + c(r, t ), 


(13) 


where 


and 


S T = -G t S — SG + SLS - K , 
c r = ~(G t -SL)c, 


S(T,t) — 0, 

c(T,t) = Hi |„ T (1 + It) - A.F. 
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to be satisfied. With this analogy, the differential equation of A (t) to be 
integrated in real time is obtained as: 



(16) 


Here, at each time t, the Euler-Lagrange equations (Eqs.(5)) are inte¬ 
grated forward along the r axis, and then Eqs.(14) are integrated backward 
with terminal conditions provided in Eq.(15). Next, the differential equa¬ 
tion of A (t) is integrated for one step along the t axis so as to estimate the 
unknown parameters from Eq.(7). If the matrix Hq q is nonsingular, the 
algorithm is executable regardless of controllability or stabilizability or the 
system. 

3.1 Convergence and Stability Analysis of NRHC Es¬ 
timation Problem 

Theorem-1: If 3 V(e, 0) such that V(e, 0) < 0 V V(0) = 0, V >0 and 
e / 0 , then the closed loop dynamics defined for the drive system x = 
Ax + f{x) + D(x)Q(t), the response system y = By + f(y) + D(y)Q(t) and 
the error system e(t) = Ae(t ) + B(f(y(t )) — f(x(t))) + D(y)Q(t) — D(x)Q(t), 
using the nonlinear receding horizon control methodology, is asymptotically 
stable (i.e. there exists a stability ball Z(0,z)) G R n V z > 0, such that for 
any given initial value, the synchronization error (e) is stable and converges 
(e —» 0) as t —> oo.) 

Proof. Let’s take into account the following Lyapunov function with respect 
to e and 0 in the form of 


t+T 

(e T Qe + Q T RQ)dr. 



(17) 


where clearly we have V(0) = 0 and V > 0 V e ^ 0. 

Here, it is possible to partition the candidate Lyapunov function as 
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where e(r) = e(r; e*(t), t) and 0(r) = 0(r;0*(t),t) are optimal values for 
system, at given time instant t. Now, consider a parameter estimator defined 
as follows: 


0(r) 


0(t; 0*(t), t) for t G [t + At, t + T], 

0 for r G (t + T, t + T + At). 


(19) 


Let e(-) = e(- : e(t + At),t + At) : [t + A t,t + T + At] 6 i?" and denote 
the corresponding error trajectory with initial condition e(t + At) = e(t + 
At;e(t),t). It is clear that 


e(r) 


e(r;e*(t),t) for r e [t + At,t + T], 

0 for r G (t + T, t + T + At). 


( 20 ) 


because e(t + T;e(t + At),t + At) = e(t + T;e*(t),t) = 0 and 0 = 0 for 
r > t + T. Since 0 is not necessarily optimal for trajectory of error e(t +At) 
from time t + At [24] , it follows that 

V(e*(t)) = - / (e T (r)Qe(r) + 0 T (r)i?0(r))dr + V(e(t + At),t + At; 0) 

1 

^ - / (e T (r)Qe(r) + 0 T (r)i?0(r))<ir + F(e(t + At)) 

( 21 ) 


so that 

1 rt+At 

V(e(t + At)) - V(e*(t)) < -- j (e T (r)Qe(r) + 0 T (r)/?0(r))dr. (22) 

Since 1/ is continuously differentiable, from the Mean Value Theorem [24] 
we have 


V(e(t + At)) — V(e*(t)) 
At 


(23) 


where p(At) G (0,1). Since 


e*(t) = e(t; e*(t), t) = e(t), 0*(t) = 0(t; e*(t), t) = 0(t) 


(24) 




and 0 is continuous at t, we have that 


lim 

At->0 + 


e(t + At) — e*{t) 
A t 


f(y(t),0(t),x(t)) = f(y*(t),Q*(t)). 


Since V e D and e are both continuous, we get 


lim 

At->0+ 


V(e(t + At)) -V(e*(t)) 
At 


V e V(e"(t))f(y’(t),6'(t),x(t)). 


(25) 

(26) 


Based on the above derivation, we also have 


ct+At. 


lim-— 

At->o+ 2 At 


< lim 


(e T {r)Qe{r) + 0 T (r)i?0(r))dr 


1 


ct+At. 


Ai—>o+ 2 At 


(e T (r)Qe(r))dr 


(27) 


= ~^e T (T)Qe( T )). 


Dividing both sides of (22) by At > 0 and taking the limit as At —* 0 + yields 
to 

\7 e V(e*(t))f(y*(t),Q*(t),y(t)) < -^(e T (r)Qe(r)) < 0. (28) 


Therefore, 


dV(e*(t)) 

dt 


= V e V r (e*(f))/(y*(f), 0*(t),x(t)) < 0. 


(29) 


Obviously, V(e(t)) < 0 and therefore V(e(t)) (namely, J(e(t))) is non¬ 
creasing. 

It is clear that J* +T (e T Qe + © 7 RQ)dr —> 0 as t —> oo, which indicates 
the asymptotic stability of the synchronization error. ■ 

With this result, we show that the non-increasing monotonicity of the 
cost function is a sufficient condition for the stability and the unknown time- 
varying parameters can be estimated through this process, without any given 
reference trajectory for the parameters. 
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4 Simulation Results 


In this section, we use two benchmark examples to verify the effectiveness of 
the proposed method. 

Example-1: We first consider the Lorenz system [37] with time-varying 
parameters.The objective is to estimate the time evolution of unknown pa¬ 
rameters using the proposed nonlinear receding horizon algorithm and show 
the applicability of proposed algorithm on a nonlinear dynamics (such as 
modified Lorenz system). 

Here, modified Lorenz system with time-varying parameters is given by 

■ - 

28xi — XiX 3 — x-2 , 

XiX 2 - |x 3 

0i{y2-yi) 

282/1 - ynj:i ~ 2/2 • 

. yiy 2 - 0 2 y 3 _ 




d x{t) 
d t 

d y(t) 
d t 


For this specihc problem, the performance index is chosen as follows: 

1 r l + T 

J = 2 J {[2/(0 ~ x ( t ')] T Q[y < y t ') - ®(^)] + Q(t') T RQ(t')}dt' , (31) 

where Q = diag(l, 0.5, 0.1) and R = diag(0.02, 0.02). 

The smooth function of time interval T, taken from [35], in the perfor¬ 
mance index is given by 


T(t) = T,(\ - e-” 1 ), 


(32) 


where Tf = 0.1 and a = 0.01. The stable matrix is chosen as A s = 50/. 
The initial states are given by 


xi(0) 


-3 


1/1(0) 


-6 

^2(0) 

= 

-3 


2/2(0) 

= 

-6 

_X 3 (0)_ 


15 


jy 3 (o)_ 


22 


(33) 


The simulation is implemented in MATLAB, where the time step on the 
t axis is 0.01s and the time step on the artificial r axis is 0.005s. Fig.l shows 
the trajectories of drive-response systems in Eqs.(30) with initial conditions 
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Figure 1: (Color online) The trajectories of states with time invariant pa¬ 
rameters where solid line denotes the trajectory of response system and dash 
line denotes the trajectory of drive system. 


given in Eqs.(33). The estimated parameters 6\ and d 2 are presented in Fig. 2, 
where Fig. 3 depicts the estimated error which clearly show the estimated 
parameters converge to their true values by using the NRHC method. 

ft is possible to see from the Figs. 1-3 that the convergence rate of the 
estimated parameter is less than 5 seconds, and is relatively fast. In a real-life 
application, this methodology (and the associated estimation results) could 
be used to provide (feedback) information regarding the unknown parameter 
after a grace period (say 8-10sec) and/or the desired level of convergence is 
achieved. 

Example-2: In the following, we consider a nonlinear system taken from 

[7]: 


d x(t) 
df 

dy(t) 

dt 


—x\ — 2xi + 6i(t) 
-2x 2 + 0 2 (t)xi ’ 

~y\ ~ 2yi + 0i (t) 
-2y 2 + 0 2 {t)yi. J ’ 


(34) 
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Figure 2: (Color online) The trajectories of estimated constant parameters 
where solid line denotes the trajectory of estimated parameters and dash line 
denotes the trajectory of their true values. 
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Figure 3: (Color online)The trajectories of estimated errors. 
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where the time-varying parameters are defined as 


Oi(t) 

02 (*) 


{ 2 + sin(f) 

2 -sin(f) + sin(2t + f) 
2 — sin(|) + sin(|) 

3 cos(0.l7rt). 


0 < f < 6vr 

6 n < t < 147T + ^ 

t > 147r + 


and the initial states are given by 


2:1(0) 


'O' 


yi(0) 


'1' 

_27 2 (0)_ 


0 

? 

»2<0) 


2 


(35) 


(36) 


For this specific problem, the performance index is chosen as follows: 


1 f t+T 

J = 2 J {[j/(0 - x ( t ')] T Q[y( t> ) - ®(0] + G(t') T RG(t')}dt' , (37) 


where Q = diag(100,110) and R = diag(0.1,0.2). The stable matrix is 
designed as A s = 10/. 

Fig.4 shows the trajectories of systems given in Eqs.(34) with initial con¬ 
ditions provided in Eqs.(36). The estimated parameters 9± and 62 are shown 
in Fig.5, where Fig. 6 depicts the estimated error. 

ft is clear from Fig.5 and Fig.6 that NRHC is able to perform as desired, 
and estimated time-varying parameters converge to their true values. 


5 Conclusions 

In this study, we considered a real-time nonlinear receding horizon based 
control algorithm to estimate unknown time-varying parameters in nonlin¬ 
ear systems. For this purpose, we demonstrated the fact there exists a Lya¬ 
punov function which could be used to achieve asymptotic stability and pro¬ 
vide stability guarantees for such an approach. Two benchmark examples 
demonstrates successful results of applicability of such concept on nonlinear 
systems. 
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Figure 4: (Color online) The trajectories of states with time invariant pa¬ 
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